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Abstract. We present a unified framework for equilibrium and nonequilibrium many-body 
perturbation theory. The most general nonequilibrium many-body theory valid for general 
initial states is based on a time-contour originally introduced by Konstantinov and Perel'. The 
various other well-known formalisms of Keldysh, Matsubara and the zero-temperature formalism 
are then derived as special cases that arise under different assumptions. We further present a 
single simple proof of Wick's theorem that is at the same time valid in all these flavors of many- 
body theory. It arises simply as a solution of the equations of the Martin-Schwinger hierarchy for 
the noninteracting many-particle Green's function with appropriate boundary conditions. We 
further discuss a generalized Wick theorem for general initial states on the Keldysh contour and 
derive how the formalisms based on the Keldysh and Konstantinov-Perel'-contours are related 
for the case of general initial states. 



1. Introduction 

In many physical situations we are interested in knowing the expectation value of some observable 
quantity of a system in or out of equilibrium. For quantum systems of many identical and 
interacting particles a very convenient mathematical object to extract this information is the 
Green's function. Let p be the density matrix which describes the system at time, say, to and 
H{t) be the Hamiltonian of the system for times t > to- The n-particle Green's function G n is 
defined according to 

G n (l . . . n; 1' . . . n') = ^ Tr [p T {^(1) . . . i> H {n)$ H {ri) . . . ^(l')}] • (1) 

In this formula 1 = (xi,t{), 2 = (a?2>£2)> etc. are collective indices for the position-spin 
coordinates x = r, a and time t, the symbol Tr denotes a trace over the Fock space, T is the 
time-ordering operator and V'-ff(j) = U(to,tj)i/)(j)U(tj,to) are field operators in the Heisenberg 
picture with respect to the Hamiltonian H (hence U is the evolution operator). The quantum 
average of a n-body operator can be calculated from the equal-time Green's function G n . 

The direct evaluation of G n from Eq. (1) is, in general, an impossible task. The first difficulty 
is brought by the Hamiltonian H = Hq + Hi nt which is typically the sum of a one-body operator 
Hq and a m-body operator H{ nt with m > 2. For H{ nt ^ the field operator ipn in the 



Heisenberg picture is a complicated object and must be approximated in some clever way. The 
second difficulty consists in taking the trace over the Fock space with a density matrix p. The 
density matrix is a self-adjoint, positive semi-definite operator with unit trace and, therefore, it 
can be written as p = e~ x / Tr[X] where X is a self-adjoint operator. For instance for systems in 
equilibrium X = (3H M with j3 the inverse temperature and H M = H — pN the grand-canonical 
Hamiltonian. To make contact with this equilibrium situation we define H M = X//3 so that 

P= e — (2) 

with Z = Tr[e — ^' ffM ]. In equilibrium Z is the partition function. In order to specify the initial 
preparation of the system we can assign either p or H M since there is a one-to-one correspondence 
between the two. If we now separate H M = Hq 1 + H^ t into the sum of a one-body operator Hq 1 
and a m-body operator H^ t with m > 2 then the trace in Eq. (1) can easily be worked out for 
-^Mnt = whereas we have to use suitable approximation schemes for H^ t ^ 0. 

Different Many-Body Perturbation Theories (MBPT) have been put forward to overcome 
these difficulties. The most popular MBPT's are probably the zero-temperature (real-time) 
Green's Function Formalism (GFF) and the finite-temperature (imaginary-time) Matsubara 
GFF [1]. These two formalisms are limited to equilibrium situations. Systems driven out of 
equilibrium by an external field are usually studied within the (adiabatic real-time) Keldysh 
GFF [2, 3]. The Keldysh GFF, however, neglects the effect of initial correlations which arc 
relevant in the short-time dynamics of general quantum systems, such as in transient dynamics 
in quantum transport or in the study of atoms and molecules in external laser fields. There 
exist two alternative GFF's to include initial correlations. The first is based on the idea of 
Konstantinov and Perel' [4] and consists in attaching the imaginary-time Matsubara track to 
the original Keldysh contour, see Refs. [3, 5, 6]. The second GFF does instead account for initial 
correlations through extra Feynman diagrams, the evaluation of which requires the knowledge 
of the reduced n-particle density matrices 

r„(xi ...x n ;x[... x' n ) = Tr^V^xi) • • • ^ {x' n )ip(x n ) • • • ip{xi)] (3) 

where the symbol Tr signifies a trace over the Fock space, see Refs. [7, 8, 9, 10, 11, 12]. These 
last two formalisms are both exact and hence equivalent. 

In all the aforementioned GFF's the dressed (interacting) G n is expanded in powers of 
the interaction Hamiltonian (H m t and/or H^ t ), leading to an expansion of G n in terms of 
the bare (noninteracting) Green's functions Go,n- The appealing feature of any GFF is the 
possibility of reducing the Gq iTI to an (anti)symmetrized product of Go = Go,i by means 
of Wick's theorem [13]. Even though the mathematical structure of all GFF's is identical, 
these formalisms are usually treated as independent probably due to the fact that the existing 
proofs of Wick's theorem are very much formalism-dependent. In this paper we show that 
Wick's theorem is the solution of a boundary problem for the Martin-Schwinger Hierarchy 
(MSH) [14] and that different GFF's correspond to different domains and parameters for the 
MSH [15]. In this way we can easily explain the common mathematical structure of every 
GFF and see how, e.g., the Keldysh GFF reduces to the zero-temperature GFF in equilibrium 
or the Konstantinov-Perel' GFF reduces to the Keldysh GFF under the adiabatic assumption. 
Our reformulation also allows us to prove a generalized Wick's theorem for interacting density 
matrices p. This naturally leads to the diagrammatic expansion with extra Feynman diagrams 
previously mentioned. The generalized Wick expansion has a form identical to that of a Laplace 
expansion for permanents/determinants (for bosons/fermions). Consequently, the calculation of 
the various prefactors is both explicit and greatly simplified. In this contribution we only state 




Figure 1. Contour 7 of Eq. (4). The contour 
consists of a forward branch going from to to 
00 (on this branch the points are denoted by 
z = t_) and a backward branch going from 00 
to to (on this branch the points are denoted 
by z = t+). 



the generalized Wick's theorem and refer the reader to Refs. [12, 15] for the proof. We will, 
however, discuss the equivalence between the GFF based on the generalized Wick's theorem and 
the Konstantinov-PereP GFF. 

2. General formula for the Green's function 

The n-particle Green's function in Eq. (1) can also be written as 



G n (l...n;l'...n') = -Tr \p T \e~ l ^ dzH{z) ^{l) . . .^{n)^ {n') . . .^') 



1' 



(4) 



Let us explain this formula and discuss the equivalence with Eq. (1). In Eq. (4) the integral 
is over the contour 7 of Fig. 1 which goes from to to 00 and back to to whereas T is the 
contour ordering operator which rearranges operators with later contour arguments to the left. 
We denote by z = t± the points on 7 lying on the lower/upper branch at a distance t from the 
origin and define the field operators with arguments on the contour as 

tp(x, z) = ip(x). (5) 

More generally every operator 0(t) with a real-time argument can be converted into an operator 
O(z) with a contour-time argument according to the rule 0(t+) = 0(t ) = 0(t). In particular 
H{t ) = H{t + ) = H(t). The reason to keep the contour argument in Eq. (4) even for operators 
that do not have an explicit time dependence (like the field operators) stems from the need of 
specifying their position along the contour, thus rendering unambiguous the action of T. Once 
the operators are ordered we can omit the time arguments if there is no time dependence. For 
instance if t\ < £2 then 

T L-i I, dzH ^^( Xutl _)^(x 2 ,t2-)} = ±U(to, oo)L>(oo, t 2 )ft(x 2 )U(t 2 ,t 1 )4>(x 1 )U(t 1 ,to) 



r{^H(xi,ti)ft(x 2 ,t 2 )}, (6) 



where the ± sign in the first equality is for bosons/fermions. One can verify that Eq. (6) is 
valid also for t\ > t 2 . This example can easily be generalized to many field operators. We 
conclude that Eq. (4) is equivalent to Eq. (1) for contour arguments on the upper branch of 
7. The G n in Eq. (4) is, however, more general since the contour arguments can lie either on 
the upper or lower branch of 7. Quantities like photoemission currents, hyper-polarizabilities 
and more generally high-order response properties require the knowledge of this more general 
Green's function. 

The density matrix in Eq. (4) can be incorporated into the contour ordering operator if we 
extend 7 as illustrated in Fig. 2 and define the Hamiltonian with imaginary-time arguments as 
H(t -i T ) = H M . Since 

e PH =e Jt W =T e J'o \ (7) 
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Figure 2. Contour 7 and Hamiltonian along the 
contour 7 to get the exact (Konstantinov-PereP) 
Green's function from Eq. (8). 



we have 



G n (l...n;l'...n') 



l Tr 


M 


e -i/ 7 dzH(z)^ _ _ _ ^( n )^t( n ') . . . ^(1') 






Tr 











(8) 



-i f H>+ dzH(z) ~ * 

where in the denominator we took into account that T{e °- } = ?7 (to, oo)U(oo, to) = 1- 

Equation (8) is, by construction, equivalent to Eq. (4). It gives the exact Green's function 
provided that the integral is done along the contour 7 of Fig. 2 and provided that the 
Hamiltonian changes along the contour as illustrated in the same figure. We now show that 
the Green's function of every GFF can be written as in Eq. (8), the only difference being the 
shape of 7 and the Hamiltonian along 7. 

We mentioned in the introduction that the calculation of the trace simplifies if the density 
matrix is of the form po = e~@ H o /Zq, with Hq 1 a one-body operator and Zq = Ti[e~^ H o ]. It 
is possible to turn a trace with p into a trace with po if the adiabatic assumption is fulfilled. 
According to the adiabatic assumption one can generate the density matrix p with Hamiltonian 
H M = Hq 1 + H^ t starting from the density matrix po with Hamiltonian Hq 1 and then switching 
on H^ t adiabatically, i.e., 



U v (t , -00) 



a ttM 

e -P H o 



U v (-oo,t ) = U v (t , -00) poU^-oc, t ), 



(9) 



Z " vu ' ' Z 

where U v is the real-time evolution operator with Hamiltonian 

H v (t) = Hf + e-^\H™, (10) 

and i] is an infinitesimally small positive constant. This Hamiltonian is equal to Hq 1 when 
t — > — 00 and is equal to the full interacting H M when t = to- in general the validity of the 
adiabatic assumption should be checked case by case. Under the adiabatic assumption we can 
rewrite Eq. (4) as (omitting the arguments of G n ) 



G n = — Tr 



Po U v (-oo, t )T |e ' ff dzH{z) 4>(l) . . . ^(n)^t( n ') . . . </,(!')} U v (t , -00) 



(11) 



We now see that if we construct the contour 7 of Fig. 3 and let the Hamiltonian change along 
the contour as 



H(t±) = < 
H{z € 7 M 



M 
int 



H(t) = H (t) + H int (t) 



for t < to 
for t > to 



F M = H - pN, 
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Figure 3. Contour 7 and Hamiltonian 
along the contour 7 to get the adiabatic 
(Keldysh) Green's function from Eq. 
(8). 



then Eq. (11) takes the same form as Eq. (8). We will refer to this way of calculating G n as the 
adiabatic formula. This is exactly the formula used by Keldysh in his original paper [2]. The 
adiabatic formula is correct only provided that the adiabatic assumption is fulfilled. 

We can derive yet another expression of G n for systems in equilibrium at zero temperature. 
In equilibrium H M = H — pN and therefore Hq 1 = Hq — pN and H^ t = Hint- Assuming that 
Ho and i?i nt commute with N the evolution operator U v in Eq. (9) can be calculated with 
Hamiltonian 

H ri (t) = H + e- r '\ t - t »\H illt (12) 

since the addition of —pN corresponds to multiplying by a phase factor. In Eq. (9) this 
phase factor cancels out since U v (— 00, to) = [U v (to, —00)]'. Furthermore, for any finite contour- 
times in G n we can approximate the evolution operator U in the field operators ipn with the 
evolution operator U v since we can always choose 77 <C l/\t — to\ and hence H v ~ H. Thus Eq. 
(11) becomes 



Tr 



4>(n)^(n')...^(l' 



(13) 



where 7 is a contour that goes from —00 to 00 and back to —00 and H(t±) = H n {t) is the 
Hamiltonian of Eq. (12). Next we observe that the interacting p can also be generated starting 
from po and then propagating backward in time from 00 to to using the same evolution operator 
U v since H v (to —t) = H v (to +t). In other words 

P = U v (to,oo) poU v (oo,to). (14) 

Comparing this equation with Eq. (9) we conclude that 

Po = U r ,(-oo,oo)p U ri (oo,-oo). (15) 

If the ground state |<J>o) of Hq — pN is nondegenerate then the zero-temperature po = |$o)( c ^o| 
is a pure state and Eq. (15) implies that 

($o|£>r,(oo,-oo) = e iao ($ \. (16) 

We will refer to the adiabatic assumption in combination with equilibrium at zero temperature 
and with the condition of no ground-state degeneracy as the zero-temperature assumption. The 
zero-temperature assumption can be used to manipulate Eq. (13) a bit more. We have 



po = |$o}($c 



l^oH^ol^^-oo) 

($ |^r,(oO,-Oo)|$o) 



lim — 

/3->00 fJTj. 



^OO, -OO ) 



e~P H o ^(00,-00) 



(17) 



Y 



Figure 4. Contour 7 and Hamiltonian 

along the contour 7 to get the zero- 
temperature Green's function from Eq. 
(8). 



00- i(3 



Inserting this result into Eq. (13) we find that the zero-temperature Green's function can again 
be written as in Eq. (8) with the contour 7 that starts at —00, goes all the way to 00 and then 
down to 00 — i(3, see Fig. 4, and with the Hamiltonian H(z) that varies along the contour as 
illustrated in the same figure. It is worth noticing that the contour 7 has the special property of 
having only a forward branch and that for the zero-temperature assumption to make sense the 
Hamiltonian of the system must be time independent. There is indeed no reason to expect that 
by switching on and off the interaction the system goes back to the same state in the presence 
of external driving fields. 

To summarize the exact (Konstantinov-PereP), adiabatic (Keldysh) and zero-temperature 
Green's functions have the same mathematical structure, given by Eq. (8). What changes is 
the contour and the Hamiltonian along the contour. 



3. Wick's theorem and Many-Body Perturbation Theory 

To be concrete we specialize the discussion to interaction Hamiltonians H mt and H^ t which are 
two-body operators. Higher order n-body operators lead to more voluminous equations but do 
not rise conceptual complications. Thus we write 



flint (*) 



^ J dxidx 2 v(x 1 ,x 2 ; z^^xi)^ (x 2 )^(x 2 )^{xi). 



(18) 



In the exact (Konstantinov-PereP) formula the interaction v(x±,x 2 ; z) = v(x±,x 2 ) is the 
interparticle interaction for z on the horizontal branches whereas v(x±, x 2 ; z) depends on the 
initial preparation for z on the vertical track. For instance in equilibrium v(x\,x 2 ;tQ — it) = 
v(x±,x 2 ). On the other hand in the adiabatic (Keldysh) and zero-temperature formula 
v(xi,x 2 ;z) = e~ v ^~ to ^v(xi,x 2 ) for z on the horizontal branches whereas v(xi,x 2 ;z) = for 
z on the vertical track. Let us consider Eq. (8) and write the exponential of H as the product 
of the exponentials of Hq and flint: 



G n (l...n; l' ...n) = 



l Tr 


v\ 


■ -i j; dzH {z) e -i f, dzH int (z)^ {1) ^ (n) ^ t („') . . . ^(1')' 




jn 


Tr 


M 


-i ( dzHo(z) —i 1 dzHi nt (z) 


[ 





(19) 



The expansion in powers of flint leads to an expansion of G n in terms of noninteracting Green's 
functions Go, n - The Go, n are obtained from Eq. (19) by setting flint (z) = for all z € 7. For 
instance for n = 1 we get 

§ F Of I !')••■ v(k\ k')G , 2k+1 (a, 1,1',...; 6, 1+, 1 /+ , . . .) 

(20) 



G(a;b) = 



k=0 



§ f Of I v ^ !') ■ ■ ■ fc')Go,2 fe (i, i', • • • ; i + , i /+ , • • •) 

k=0 ' v 7 



for n = 2 we get 

E f (I) / !') ■ ■ ■ «fo fc , )G 0>2 fc+2(a, 6, 1, 1' . . . ; c, d, 1+, l'+, . . .) 
G 2 (a,6;c,d) = — ~ k , (21) 

E F (!) / !') • • • fc , )G 0) 2fc(i, i+, • • •) 

etc. In these equations a = (x a , t a ), b = (xb, tj,) are collective indices like 1, 2, . . ., the interaction 
v (j> f) = ^(^ii z 'j) v ( x i-> x 'j'i z j) an d the integrals are over 1,1',..., k, k' . 

The appealing feature of any GFF is the possibility of reducing the noninteracting Go, n to 
a (anti)symmetrized product for (fermions) bosons of one-particle Green's functions Go- This 
reduction is called Wick's theorem. The existing proofs of Wick's theorem are rather laborious 
and differ depending on whether one is working with the zero-temperature or Matsubara or 
Keldysh Green's functions. Below we give a simple and general proof of Wick's theorem which 
applies to all cases. 

We consider a one-body Hamiltonian of the form 



Hq(z) = J dx $(x)h(x, z)ip(x). 



(22) 



The more general case of a nondiagonal h(x,x',z) can be treated in a similar manner. The 
Green's functions Go, n satisfy the noninteracting MSH 



i^- - h{k) 
dzk 



G 0) „(l ...n;l'...n') = E (±) k+j S(k; f) G 0) „_i(l ... fc ... n; l' n') 

i=i 



G ,„(l...n; l'...n') 



(23) 



J] fc') Go,n-i(l ... j ... n; 1' n') 

(24) 

where the hook over the arguments in Go, n -i means that those variables are missing. The MSH 
is a set of coupled differential equations to be solved on the contour 7 of the Green's function 
of interest (exact, adiabatic, or zero-temperature). In all cases from the definition Eq. (8) 
it follows that the Gq^ satisfy the Kubo-Martin-Schwinger (KMS) relations, i.e., the Go, n are 
(anti)periodic along the contour 7 with respect to all their contour arguments. Therefore we 
can calculate the Go, n by solving the MSH with KMS relations. We now show that the solution 
is given by the Wick theorem 



G ,„(l, • • • ,n;l',... ,ri) 



Go(l;l') 



Go(l;n') 



G (n;l') ... G (n;n') 



(25) 



± 



where the symbol | . . . \± signifies the permanent/determinant for the case of bosons/fermions 
and Go is the solution of Eqs. (23) and (24) with n = 1, i.e., 



G (l;l') = 5(l;l / ), G (l;l') 



5(1;!') 



(26) 



with KMS boundary conditions. Expanding the permanent/determinant along row, say, k we 
get 

n n n 

G ,„(l, ...,n;l',...,n') = ]T(±) fc ^G (£:, /)G 0) „-i(l ... jfe ... n; l' n') (27) 

i=i 



which is clearly a solution of Eq. (23). Similarly, we can readily verify that Eq. (25) is also 
solution of Eq. (24) by expanding the permanent/determinant along column k. It remains to 
check that the Go jra in Eq. (25) fulfills the KMS relations. The contour argument z/, appears in all 
the Go of the k-th row of Eq. (25) and nowhere else. Therefore when we move from the starting 
to the ending point of 7 all entries of row k pick up a (±) sign. Since the permanent/determinant 
of a matrix in which we multiply a row by (±) is (±) the permanent/determinant of the original 
matrix we conclude that Go, n is (anti)periodic with respect to the first n contour arguments. 
With a similar reasoning one can prove that Go, n is (anti)periodic with respect to the last n 
contour arguments. This concludes the proof. 

The Wick theorem has been proven without any assumption on the shape of the contour 
and without any assumption on the form of the single-particle Hamiltonian h(x, z) along the 
contour. Inserting Eq. (25) into, e.g., Eq. (20) we get the MBPT formula for the one-particle 
Green's function 



G(a; b) 



£ h(i) k fv(l;l')..v(k;V) 



k=0 



G (a;6) G (a;l+) ... 
G (l;6) G (l;l+) ... 



G (k';b) G (fc';l +> 



g £ v)..v{k-k') 



k=0 



G (l;l+) Go(l;l'+) 
Go(l';l+) G (T;1'+) 



G (a;fc'+) 
Go(i;* ,+ ) 

G (k';k' + ) 



G (l;k'+) 
Go(l';fc /+ ) 



G (k';l+) G (k'; 1'+) ... G (k';k'+) 



± 



(28) 



which is an exact expansion of the interacting G in terms of the noninteracting Go- The MBPT 
for higher order Green's functions can be derived similarly. In the next Section we discuss how 
the variuos GFF's follow from Eq. (28). 



4. Matsubara, Keldysh and zero-temperature formalisms 

In the Konstantinov-PereP formalism the Green's function G is given by Eq. (28) where the 
z- integrals run on the contour of Fig. 2. It is worth stressing that if the times t a and tf, in 
G(a; b) are smaller than a maximum time T then it is sufficient to perform the z integrals over 
a shrunken contour like the one illustrated in Fig. 5. This is a direct consequence of the fact 
that if the contour is longer than T then the terms with integrals after T cancel off [15]. 




Figure 5. The shrunken contour 7 that 
can be used in Eq. (28) to obtain the 
Green's function with real-time smaller 
than T. 
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The Matsubara GFF is used to calculate Green's function with imaginary times and it is 
typically applied to systems in equilibrium at finite temperature. For this reason the Matsubara 



GFF is also called the "finite-temperature formalism". To calculate G with imaginary-time 
arguments we can choose T = to in Fig. 5 and hence shrink the horizontal branches to a point 
leaving only the vertical track. Therefore the Matsubara GFF consists of expanding the Green's 
function as in Eq. (28) with the z-integrals restricted to the vertical track. It is important to 
realize that no assumptions, like the adiabatic or the zero-temperature assumption, are made 
in this formalism. The Matsubara GFF is exact but limited to initial (or equilibrium) averages. 
Equivalently we can say that the Matsubara G is the same as the Konstantinov-PereP G on the 
vertical track. 

The formalism originally used by Keldysh was based on the adiabatic assumption. The 
Keldysh Green's functions are again given by Eq. (28) but the z-integrals are done over the 
contour of Fig. 3 and the Hamiltonian changes along the contour as illustrated in the same 
figure. The important simplification of the Keldysh GFF is that the interaction v is zero on 
the vertical track. Consequently in Eq. (28) we can restrict the z-integrals to the horizontal 
branches. Like the Konstantinov-PereP formalism, the Keldysh GFF can be used to deal with 
nonequilibrium situations in which the external perturbing fields are switched on after time to. 
In the special case of no external fields we can calculate interacting equilibrium Green's functions 
at any finite temperature with real-time arguments. 

The zero-temperature formalism relies on the zero-temperature assumption. As we already 
discussed this assumption makes sense only in the absence of external fields. The corresponding 
zero-temperature Green's function is given by Eq. (28) in which the z-integrals are done over the 
contour of Fig. 4 and the Hamiltonian changes along the contour as illustrated in the same figure. 
Like in the Keldysh GFF the interaction v vanishes along the vertical track and hence the z- 
integrals can be restricted to a contour that goes from — oo to oo. The contour ordering operator 
is then the same as the standard time-ordering operator. For this reason the zero-temperature 
Green's function is also called time-ordered Green's function. The zero-temperature GFF allows 
us to calculate the interacting G in equilibrium at zero temperature with real-time arguments. 
It cannot, however, be used to study systems out of equilibrium and/or at finite temperature. In 
some cases, however, the zero-temperature formalism is used also at finite temperatures (finite 
(3) as the finite temperature corrections are small. This approximated formalism is sometimes 
referred to as the real-time finite temperature formalism [1]. We emphasize that in the real-time 
finite-temperature formalism (like in the Keldysh formalism) the temperature enters in Eq. (28) 
only through Go which satisfies the KMS relations. In the Konstantinov-PereP formalism, on 
the other hand, the temperature enters through Go and through the contour integrals since the 
interaction is nonvanishing along the vertical track. 

5. Generalized Wick's theorem 

The MBPT of the exact GFF requires the knowledge of the operator H M on the vertical track. 
In many physical situations, however, it is easier to specify the initial state (or initial density 
matrix) instead of H M . In these cases the preliminary step to apply MBPT consists in obtaining 
H M from p, something that can be rather awkward. For instance if p = \ is a pure state 

then H M is an operator with |^) as the ground state. The existence of a generalized Wick's 
theorem that uses directly p would be of very valuable. In this Section we will show how to 
construct such a generalized framework. 

We consider again Eq. (4) but this time we do not incorporate p in the contour ordering. 
In Eq. (4) the contour 7 is that of Fig. 1 and the Hamiltonian is the physical Hamiltonian 
which, for simplicity, we take as the sum of Ho in Eq. (22) and Hi nt in Eq. (18). We write 
the exponential in Eq. (4) as the product of two exponentials, one containing Ho and the other 
containing H{ nt , like we did in Eq. (19). The subsequent expansion of G in powers of H mt leads 



to the expansion 



00 1 /i\ k r 



(29) 



and similarly for higher order Green's function. In Eq. (29) the Green's functions g n are 
noninteracting Green's functions averaged with an arbitrary density matrix p 



g n (l...n; l' 



1 



. n 



Tr 



pT{e * A ^^(l) . . . $(n)ft(ri) . . . ^(1'; 



(30) 



We will now prove a generalized Wick theorem to write these g n in terms of the one-particle 
Green's function g = g\ and the n-particle reduced density matrices T n defined in Eq. (3). 

The Green's functions g n satisfy the noninteracting MSH on the contour of Fig. 1. The 
problem in solving the MSH to obtain the </ n 's is that we cannot use the KMS relations as 
boundary conditions. Indeed it is easy to verify that the g n are not (anti)periodic along the 
contour. A convenient choice of boundary conditions follows directly from the definition of g n 
and reads 

{±i) n lim g n (l . . . n; l', . . .ri) = T n (x 1 . . .x n ;x[ . . .x' n ) (31) 

Zk,Zj->to- 

where the limit is taken with the order z\ < . . . < z n < z' n < . . . < z[ of the contour arguments. 
The permanent/determinant 



g n (l...n;l' ...ri) 



0(1; 1') 

g(n; V) 



g(l;n') 
g(n;n') 



= |5| n (l...n; I' ...ri) 



(32) 



± 



is a solution of the MSH but, in general, with the wrong boundary conditions. In Eq. (32) the 
symbol | . . . |± signifies the permanent/determinant of the matrix inside the vertical bars. The 
particular solution must be supplied with the solution of the homogeneous equations 



iJ- - h{k) 
dz k 



n (l...n;l'...n') = O 



(33) 



g n (l...n;\' . 



. n 







(34) 



to satisfy the correct boundary conditions. We observe that g n is not discontinuous when its 
contour arguments cross each other since in the right hand side of Eqs. (33) and (34) there is no 
(5-function. Consequently the equal-time limit of g n is independent of the order of the contour 
arguments. 

Let us start by showing how to solve the MSH for g^. We write 52 = 1 5 1 2 + §2 where g satisfies 
the first equation of the MSH with boundary conditions 



(±i) lim g(l-l') = r 1 (x 1 ;x' 1 ) = r(x 1 ;x' l ). 

The boundary conditions for g 2 follow directly from Eq. (31) and read 

(±if lim g 2 = (±if lim (g 2 - \g\ 2 ) = T 2 - |r| 2 = C 2 . 



(35) 



(36) 



Next we consider the spectral function on the contour 

A(l;l') = i[g>(l-l')-g<(l;l')]. (37) 
This function takes the same value for z\ = t\± and z[ = t[± and satisfies the equations 



.4(l;l') = ^(l;l') 



= 0. 



(38) 



Furthermore, due to the (anti) commutation rules of the field operators 

A(xi, z; x[, z) = 5(xi — x[). (39) 

Therefore 

g 2 (l, 2; l', 2') = ^ dx\dx 2 dx' x dx' 2 A(\\ x\, t -)A(2; x 2 , to-)C 2 (xi, x 2 ; x[, x' 2 ) 

x A(xi,t O -;l0^2,io-;2') (40) 

is clearly the solution of the homogeneous MSH with the correct boundary conditions, see Eq. 
(36). We can manipulate Eq. (40) by introducing a linear combination of 5- functions on the 
contour 

6-(z)=S(z,t -)-S(z,t 0+ ). (41) 
The spectral function appearing in Eq. (40) can be written as 

A(l;x' 1 ,tQ-) = i j dz g(l; x[, z)5-(z) 
J-f 

and 



(42) 
(43) 



A(x 1 ,t ^) = - i [d,6^g(x 1 ,^). 
Inserting these expressions into Eq. (40) we find 

52(1, 2; 1', 2') = J dld2dVd2'g{l- l)g(2; 2)C 2 (1, 2; V, 2')g(l>; l')g(2'; 2') (44) 

where 



C 2 (l,2;l',2>) = 5^z 1 )5-(z 2 )C 2 (x 1 ,x 2] x[,x' 2 )5„(z[)5-(z' 2 ) 



(45) 



is the two-particle initial-correlation function. In conclusion g 2 can be written in terms of g and 
T„ with n < 2. Since g 2 = \g\ 2 + g 2 this result provides a decomposition of g 2 in terms of g and 
reduced n-particle denity matrices. 

The generalization of Wick's theorem to g n reads 



n-2 



9n = \g\n+ ££(±) |P+QI \g\i(P;Q')~gn-i(P;Q') + ~g n . 

1=1 PQ 



(46) 



In this formula P and Q are a subset of I ordered indices between 1 and n whereas P and Q 
is the ordered complementary subset. For instance if n = 3 and I = 1 then we can have P = 1 
and hence P = (2, 3), or P = 2 and hence P = (1, 3), or P = 3 and hence P = (1, 2). The sign 
of the various terms is given by |P + Q\ = J2i=i(Pi + Qi) where pi and % are the indices in the 



/-tuple P and Q. The solution of the homogeneous MSH with the correct boundary conditions 
is 

~g n {\ . . . n; l' . . . ri) = J g(l; 1) . . . g(n; n)C n (l . . . n; l' . . . n')g(n'; n') . . . g(V; l') (47) 

where the integral is over all barred variables and the n-particle initial-correlation functions are 
given by 

C n {\ ...n-l' ...n') = 8.{ Zl ) . . . 5^{z n )C n {x l . . . x n ; x\ . . . x' n )5-(z[) . . . 5„{z' n ) (48) 

with 

n-2 

C n = T n - £ Y,( ± ) lP+Ql \ T \i(Xp;X' Q )C^ l (X p -X^ - \T\ n . (49) 
1=1 PQ 

This is a recursive formula for the C n . The collective coordinate Xp = ( subset of 

the coordinates (x\ . . . x n ) and similarly the collective coordinate X'q = (x' qi . . . x' qi ) is a subset 
of the coordinates (x[ . . . x' n ). 

We defer the reader to Ref. [12] for the proof of the generalized Wick theorem. Here we 
observe that with the generalized Wick theorem we can express g n in terms of g and Y n . Since 
T n can easily be calculated from p the generalized Wick theorem is especially suited to do MBPT 
when we know p instead of H M . We further observe that the generalized Wick theorem has the 
same mathematical structure of the Laplace expansion for the permanent/determinant of the 
sum of two matrices A and B [15] 

m—l 

\A + B\ m = \A\ m +Y J T.^ P+Q ^P-,Q)\B\ m -i{P-,Q) + \B\ m (50) 
1=1 PQ 

where |^4|;(P;Q) is the permanent/determinant of the / x I matrix obtained with the rows P 
and the columns Q of the matrix A. The same notation has been used for the matrix B. With 
the identification A^j = g(k;f) and \B\ m _i(P;Q) = g m ~i(P;Q') for I = 1 . . . m — 2 and the 
definition g\ = Eqs. (46) and (50) become identical. We can thus symbolically write the 
generalized Wick theorem as 

9m = \9 + g\m (51) 
whose precise meaning is given by Eq. (46). 

6. Relation with the Konstantinov-Perel' formalism 

In Eq.(29) we have seen how we can expand the Green's function G into noninteracting Green's 
functions g n satisfying the generalized Wick theorem (46). This can be used to define a 
diagrammatic expansion of the Green's function. Let us see what terms we get when we insert 
Eq.(46) into Eq.(29). The first term \g\ n in Eq.(46) generates the usual series of connected 
diagrams for the Green's function in powers of the interaction (the disconnected diagrams are 
zero since they are integrated from to to to an d have no external points). The remaining terms 
in Eq.(46) are linear in the functions g m with m = 2, . . . ,n. As a consequence of Eq.(47) these 
functions can be diagrammatically represented by blocks C m with m ingoing and m outgoing 
g lines. Since the Wick expression (46) is linear in the g m each diagram contribution to the 
Green's function contains at most one correlation block. In Fig. 6 we display the diagrams for 
the Green's function up to first order in the interaction involving four diagrams with a C2 block 
and one diagram with a C3 block. The last three diagrams, however, vanish since for those 
diagrams there are internal time-integrations for Green's funtion lines that enter and leave a C m 
block that can be reduced to a point, since the initial correlation block only exists at time to-. 
The fact that the C m blocks are not repeated within a single diagram prevents us from deriving 




+ " + 2 



Figure 6. Diagrammatic expansion of 
the Green's function to first order in the 
interaction. 



an irreducible self-energy in terms of them. We can, however, define a reducible self-energy ay 
and write the Green's function as 

G(l;2)=0(l;2)+ / dld2g{l; i> r (l; 2)g{2; 2) (52) 

where we used the notation 7 for the contour of Fig.l to distinguish it from the extended contour 
7 that we will use later. The reducible self-energy can be split into three contributions. These 
are the sets of self-energy diagrams that start with an interaction line at their entrance vertices, 
and end with a correlation block at their exit vertices, denoted by 07", the diagrams that start 
with a correlation block and end with an interaction line, denoted by of, and the remaining 
diagrams (which either contain no correlation block or contain a correlation block only attached 
to internal vertices) which we will denoted by ay. The labels L and R therefore refer to the 
location of an interaction line at the entrance vertices. We can thus write 

a r = a r + of + of (53) 

where 

a r L (l;2) = a^(l;x 2 )5.(z 2 ) (54) 

a r R (l;2) = (5_(ziK R (£ i; 2). (55) 

We emphasize that ay is a reducible self-energy and therefore it should not be confused with 
the self-energy of Ref. [7]. We will now proceed to connect the reducible self-energy oy to 
the irreducible self-energy appearing in the equation of motion for the Green's function on 
the extended contour. As discussed in the introduction the initial ensemble is of the form 
p = e~ x /Tr[e^ x ], where X = Y^ m X m is in general a sum of m-body operators of the form 

X m = — [ dxi . . . dx' m v m (xi, . . . , x m , x[, . . . , x' m )^ (xi) . . . ^ (x m )^(x' m ) . . . i!(x[) (56) 
nr. J 

The functions v m must be chosen in such a way that the pre-scribed density matrices of Eq.(3) 
are obtained. This is, in general, a difficult task. If we, however, assume that we have succeeded 
in this task we can define the Matsubara Hamiltonian H M = X//3 and use the contour of Fig. (2) 
and the expression for the Green's function (8) to expand the Green's function diagrammatically 
in powers of Go using the standard Wick theorem of Eq.(25) [3, 5, 6]. The diagrammatic rules 
for the Green's function in the case of m-body operators are given in Ref. [5]. If we collect 
the irreducible pieces of this expansion in an irreducible self-energy S then we obtain a Dyson 
equation on the extended contour 

G(l; 2) = G (l; 2) + / d3dA G (l; 3)S(3; 4)G(4, 2). (57) 

It will be more convenient to write this as equations of motion 

(id Zl -fc(l))G(l;2) = <5(1;2)+ /d3S(l;3)G(3;2) (58) 

G(l;2)(-iK 2 -h(2)) = 5(1; 2) + /d3G(l;3)S(3;2) (59) 



When split into various components these are simply the Kadanoff-Baym equations on the 
extended contour. Let us now write the contour 7 as 7 = 7©7 M , where 7 M denotes the vertical 
(or Matsubara) track of the contour and 7 the remaining piece, which is identical to the contour 
of Fig. 1. The last term on the r.h.s of Eq.(58) can then be written as 

/ d3E(l;3)G(3;2) = / d3 E(l; 3)G(3; 2) - i [ dx 3 [ dr S 1 (1; x 3 t)G^(x 3 t; 2) (60) 

J-y J7 J JO 

where we introduced the parametrization z = to — ir on the vertical track 7 M . For a general 
function A(z, z') on the contour (spatial coordinates suppressed) we further defined 

A\z = t ± ,r) = A(z = t±,t - ir) (61) 
A^(r,z = t±) = A(t -iT,z = t±). (62) 

From the Dyson equation (57) and the Langreth rules on the contour 7 we can further derive 
that [15, 16] 

G r (l; 2) = -ij dx G M (1; xt )G A (xt ; 2) + [G M * E r • G A ](1; 2) (63) 

where * denotes a convolution between to an d to — i(3 on the vertical track and • denotes 
a convolution between to an d 00. We further defined the advanced and Matsubara Green's 
functions as 

G A (1;2) = -0(t 2 -t 1 )[G>(l;2)-G<(l;2)] (64) 

G M (1;2) = G(x 1 t -iT 1 ;x 2 to-iT 2 ) (65) 

If we further use that 

[ dz5-(z)G(xz;xt±) = G A (xt ;xt) (66) 
we find by inserting Eq.(63) into the last term of Eq.(60) that 

-i J dx 3 d T ^(l;x 3 T)G^x 3 T;2) = 

/(f3[El *G M *Er](l;3)G'(3;2) -i f d3 [E^ ★ G M ](1; x 3 t )5-(z 3 )G(3; 2) (67) 
If we therefore define El by 

E L (1;2) = -i[sUG M ](l;x 2 to)L(z 2 ), (68) 
we can rewrite the equation of motion (58) for the Green's function as 

(id zi -/i(l))G(l;2) = 5(1; 2) + / d3 [E + EU G m ★ E^ + E L ](1; 3)G(3; 2). (69) 
A similar procedure can be carried out for the adjoint equation (59). We find 

G(l;2)(-^ -h{2)) = 5(1; 2) + / d3 G(l; 3)[E + EU G M ★ E^ + S R ](3; 2) (70) 
where we defined 

E R (1;2) = -iS-iz^^idl^to^). (71) 



Given the self-energy and the Green's function with arguments on the imaginary track 7 M , we 
can regard Eqs. (69) and (70) as equations of motion for the Green's function G on the contour 
7. These equations can be integrated using the noninteracting Green's function g of Eq.(35) 
since it satisfies 

(±i) lim g(l;l') = (±i) lim G(l; l') = V{ Xl - x[). (72) 

Zl,z[^>-to- Z\,z\-+t§- 

If we, therefore, define the total self-energy as 

S tot = E + sUG M *S^ + S L + E R (73) 
we can write G in terms of two equivalent Dyson equations 

G(l;2) = 2)+ /_d3d4 5 (l;3)I! t ot(3;4)G(4;2) (74) 

G(l;2) = 5 (1; 2)+ /_ d3d4 G(l; 3)S tot (3; 4) 5 (4; 2) (75) 

To check that these equations are equivalent to the Eqs. (69) and (70) we need to be careful. 
The standard approach is to act with the operator of the form id z — h and its adjoint on both 
Dyson equations and use the equation of motion for g 

(id zi -h(l))g(l;2) = 6(1; 2) (76) 
g(l;2)(-iX 2 -h(2)) = 6(1; 2) (77) 

We need to be careful, however, since we cannot change integration and differentiation in the 
presence of delta-functions under the integral sign. The relevant integrals over the delta-functions 
need to be done first before we use Eqs. (76) and (77). In Eq.(74) we have an integral of the 
form 

jTd3<7(l;3)£ R (3;4) = -i J dx 3 [g>(l; x 3 t ) - g < (l; x 3 t )}[G M * ^}(x 3 t ; 4) (78) 

On the right hand side of this equation we recognize the contour spectral function of Eq.(37) 
that satisfies Eq. (38). We therefore see that 

(id Zl -h(l)) / d3<7(l;3)£ R (3;4) =0 (79) 

Similarly we have 

Qf d3 S L (3; 4)5(4; 2)) (-%K 2 - h(2)) = (80) 

Then by acting with id Zl — h(l) on Eq.(74) we see that we recover Eq.(69). Similarly by acting 
with —u)z~ 2 — h(2) from the left on Eq.(75) we recover Eq.(70). It only remains to check that the 
Dyson Eqs. (74) and (75) satisfy the correct boundary conditions. Since in the limit zi, Z2 — > to- 
the contribution for the integrals on the r.h.s. of the equations vanish we see that the condition 
(72) is indeed satisfied. 

Now we are ready to discuss the connection between the formulation based on the initial 
correlation blocks and the formalism based on integrations along the imaginary track. By 
comparing Eq.(74) to Eq.(52) we see that 

Cr = Stot + £tot5^tot + Stot5£tot5£tot + • • • = Stot --, ^ — (81) 

1 - 5^tot 



and hence 

Stot = 0>(1 - fl^tot) = cr r - <7 r gcr r + a r ga r ga r — . . . = a r — (82) 

l + ga r 

This yields the expansion of the irreducible self-energy Etot m terms of the Green's functions 
g and the correlation blocks C m . We would like to mention that S tot could in principle be 
calculated from the appropriate extention of the Hedin equations to include initial correlations 
[8]. This would lead to an expansion of S tot in terms of the dressed Green's function G and 
correlation blocks. However, the iterative solution of these equations depend on the starting 
point. In particular if we start with a self-energy which contains only a C2-block then the 
iterative procedure cannot generate diagrams with C-blocks of higher order. 

7. Conclusions 

We presented a unified framework for equilibrium and nonequilibrium many-body perturbation 
theory. The most general formalism for nonequilibrium many-body theory for general initial 
states is based on the Keldysh contour to which we attach a vertical track describing a general 
initial state. This idea goes back to the works of Konstantinov and PereP (who considered 
equilibrium initial states), Danielewicz and Wagner. On this contour we can straightforwardly 
prove a Wick theorem by solving the noninteracting Martin-Schwinger hierarchy for the 
noninteracting many-body Green's functions with KMS boundary conditions. This short proof 
of Wick's theorem does not need any of the usually introduced theoretical concepts such as 
normal ordering and contractions. The statement is simply that the noninteracting m-particle 
Green's function is a determinant or permanent of one-particle Green's functions. We showed 
how the various other well-known formalisms of Keldysh, Matsubara and the zero-temperature 
formalism can be derived as special cases that arise under different assumptions. We further 
discussed a generalized Wick theorem for general initial states on the Keldysh contour. It again 
arises as a solution of the noninteracting Martin-Schwinger hierarchy for the noninteracting 
many-body Green's functions but this time with initial conditions specified by initial m-body 
density matrices. The final result of Eq.(51) is an elegant alternative to the Wick theorem of 
Eq.(25) for KMS boundary conditions. We finally showed how the formalisms based on the 
Keldysh and Konstantinov- Perel'-contours are related for the case of general initial states. 
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